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INTRODUCTION 


1.0 Background 

For years much of the fracture mechanics research, both theoretical 
and experimental, has been devoted to the study of problems in which 
the loads were applied in a direction perpendicular to the plane of the 
crack. Studies of this type provide valid information for solely mode 
one, or opening mode problems. 

A knowledge of only mode one displacements is inadequate or 
insufficient to perform the proper design analyses required for the 
structures being considered today. Factors such as complex loadings, 
nonuniform loadings, and the nonalignment of loadings and crack geome- 
tries demand that a study of mixed mode fracture mechanics be accom- 
plished to insure that the maximum possible structural integrity is 
attained. Mixed mode problems arise frequently in bridges, pipes, 
pressure vessels, and aircraft and aerospace structures. The necessity 
for using high strength alloys in the fabrication of many structures 
of this type intensifies the requirement for adequate mixed mode 
fracture mechanics analysis since many of these alloys are often quite 
brittle and sensitive to the presence of flaws or mechanical defects. 

An example of a situation in which fracture mechanics considera- 
tions will be exercised is in the design of the space shuttle structure. 
Many components of the space shuttle will be reused a number of times 
and be exposed to hostile environmental effects, so a proper design 
analysis must be made which will include adequate assurance that the 
desired life is not threatened by possible crack growth. Mixed mode 
fracture mechanics considerations would be of less consequence if the 
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i 

i 

structural components were to be used only once, but they gain : 

increased magnitudes of importance when the components must withstand 
a number of launches and recoveries. It is for considerations such as 

I 

these that a study of mixed mode fracture mechanics problems was | 

i 

undertaken. 

Part I of this report is devoted to the study of a semielliptical 
surface flaw which has nonuniform shear loading applied to the surfaces 
in either of two directions, parallel or perpendicular to the major axis 

1 

of the ellipse. An alternating method is used to calculate the mode two j 

and mode three stress intensity factors. These stress intensity factors j 

are presented as a function of position along the crack border for a | 

number of crack shapes and crack depths. No similar results have i 

i 

heretofore been obtained for this problem. j 

A number of studies have been made for an elliptical crack in an 
infinite solid and a semielliptical surface flaw in a finite solid, both 
of which are subjected to uniaxial tension. However, no valid solution 
exists for a semielliptical surface flaw having small depth to length 
ratios, a/2c, which penetrates the solid with large depth to thickness 
ratios, a/t. Kobayashi (1)* has stated that semielliptical surface 
flaws with a/2c = 0.12 to 0.2 with depth ratios a/t = 0.9 are probably 
the most critical problems in applied fracture mechanics without an 
adequate solution. Part II of this report is devoted to the study of 
semielliptical surface flaws subjected to normal loading. Results of 
this study were presented in Technical Report No. 4 but are included 


* 


Numbers in parentheses refer to references at the end of the report . 
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here to allow the presentation of additional computations made since 
the publication of that report. 

2.0 Previous Work 

Inglis (2) in 1913 presented the initial analytical work leading 
to the present day concepts of fracture mechanics. He studied the 
stresses in the vicinity of an elliptical hole which, in the limit, 
becomes a sharp crack. This work then provided Griffith (3) with the 
background recessary to postulate his theory of crack propagation based 
on a balance between released elastic strain energy and absorbed 
surface energy. 

Much of the initial fracture mechanics work was done on two- 
dimensional problems, but in 1945 Sneddon (4) formulated the axi- 
symmetric three-dimensional problem of a circular flat crack in an 
infinite solid with uniform normal pressure applied to its surfaces. 

In 1960 Green and Sneddon (5) solved a similar problem for an ellipti- 
cal crack embedded in an infinite solid, and in 1971 Shah and Kobayashi 
(6) extended this solution to the case of an embedded elliptical crack 
under arbitrary normal loading. All of these solutions dealt only 
with embedded cracks. 

The solution of Green and Sneddon (5) was used by Irwin (7) in 
1966 to estimate stress intensity factors for a semielliptical crack in 
the surface of a flat plate. Smith (8, 9) in 1966 presented the 
solution to the problem involving a semicircular surface crack located 
at a frei surface of a half-space and loaded with an arbitrary normal 
pressure. He used an alternating technique originally used by 
Lachenbruch (10) to remove and calculate the effects of normal stresses 


on the free surface. 
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Smith ami Alavi (ll, 12, 13) later solved the problem of a 
circular crack oii.Soddcd in a semi - infinite solid, and , . oblea of a 
part-circular crack in the surface of a semi-infini solid. "hresher 
and Smith (14, 15) in 1972 extendeJ this work to solve the problem of a 
part-circular surface flaw in z finite thickness solid. In 1973, Shah 
and Kobsyashi gave solutions to the problems of an embedded elliptical 
crack approaching the free surface of a semi-infinite solid subjected 
to uniform tension (16, 17) anu linearly varying pressure (18). All of 
the solutions mentioned above were for normal loading of the crack 
surface, and thus involved only mode one displacements and stress 
intensity factors. 

Much less work has been done for the three-dimensional problems 
involving cracks which have shear loading applied to their surfaces, 
although some work has been done for cracks embedded in infinite 
bodies. Segedin (19) in 1950 presented a solution for an embedded 
circular crock subjected to a uniform shear stress. This sane problem 
was also solved by Westmann (20), while r.:-helby (21) discussed it as a 
special case of his more general solution for an ellipsoidal inclusion 
in an infinite elastic medium. 

Kassir and Sih (22) in 1966 formulated the problem of an embedded 
elliptical crack in an infinite solid, subjected to uniform shear 
applied in any direction. They studied the stress field near the crack 
border, and presented mode two and mode three stress intensity factors 
as functions of position along the crack border. 
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3 . 0 Objectives of the Report 

The objectives of this report are two-fold. The first objective 
is to present the stress intensity factors for a semielliptical surface 
crack in a finite thickness solid subjected to shear loading. 

The second objective is to present stress intensity factors for a 
semielliptical surface crack which has normal loading applied to its 
surfaces. This work is the first presentation of an analysi* of this 
problem which includes all known elastic effects, and it covers ranges 
of crack depth and shape for which results of direct analysis arc not 
available. 
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PART I --SHEAR LOADING OF A SEMIELLIPTICAL SURFACE CRACK 

Chapter 1 

THE SHEAR LOADING PROBLEM 

1 . 0 General Stat e n ent of the Surface Crack Problem 

Consider a semielliptical surface crack located in the front 
surface of a finite thickness slab, as shown in Fit, 1 . The crack is 
considered to lie on the plane defined by z = 0, and has a semiminor 
axis "b" located along the y-axis, and a seraimajor axis "a" located 
along the x-axis. 

The problem is to determine the stress intensity factor along the 

periphery of the crack when the surface of the crack is subjected to 

nonuniform shear stresses t and t , the positive directions of which 

zx zy K 

are shown in Fig. 1. Boundary conditions for this problem are deter- 
mined from a consideration of the conditions of skew symmetry of the 
displacements of tht problem. These boundary conditions are: 

o zz =0, (z * 0) (1.1) 



where u and v are the displacements in the x and y directions, 
respectively. The front and back surfaces of the boay are considered 
to be free of stresses. 

In addition, the shear stresses on the surface of the crack can be 
expressed by the following series representations: 
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3 3 

t = V V A x m y n , 

zx /j Z_j mn 7 ’ 

m=0 n=0 

^ < 1 ; z = o\ ( 1 . 3 ) 

\y-Zt B »n • 

m=0 n=0 

where m + n <_ 3 . The sum, m + n, is limited to being less than or 
equal to three because of the enormous quantity of work necessary to 
derive the equations for the stresses and stress intensity factors. 

2.0 The Method of Solution 

The method of solution which will be used is an alternating 
technique similar to the one which has been used by Thresher and Smith 
(15), and in a slightly different manner by Shah and Kobayashi (17). 

Two elastic solutions are required to apply this technique. The 
results from each solution are superposed to yield a solution to the 
problem. These two solutions will be referred to as Solution 1 and 
Solution 2 as follows: 

Solution 1. The determination of the stresses at any point in 
an infinite solid due to the presence of an elliptical 
crack in the infinite solid. The crack surface is subjected 
to variable shear stresses. 

Solution 2. The determination of the stresses within a semi- 
infinite body subjected to uniform normal and shear stress 
applied over a rectangular portion of the surface. 
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2. 1 Elastic Solution of a Crack in an Infinite Solid 

Solution 1 represents the elastic solution of an elliptical crack 
in an infinite solid subjected to nonuniform shear stresses. These 
shear stresses can be applied to the surface of the crack parallel to 
the minor axis or the major axis, or parallel to both axes simulta- 
neously . This method solves the three-dimensional equations of 
elasticity for the displacements and the stresses at any particular 
point within the infinite body. Additional derivations are made for 
the equations which define the stress intensity factors as functions of 
position along the crack border. The derivations and the explanation 
of this method are accomplished in Chapter 2. 


2 . 2 The Half-Space Solution 

Solution 2 utilizes the general formulation of a half-space 
problem by Love (23, 24). This is a solution for stresses in a half- 
space due to normal and shear tractions applied to a rectangular 
portion of the surface. 

For the problem under consideration here, the shear stresses t 
and t had to be calculated from the Love formulation on the face of 

zy 

the crack. This had not been done before so it was necessary to derive 
the equations using Love's general formulation as a beginning point. 

The skew symmetry conditions of this problem preclude the presence of 
any residual normal stresses o zz on the face of the crack, so it was 
not necessary to calculate this component. 

The details of the solutior are presented in Chapter 3. 


2.3 The Alternating Method 

The alternating method which is applied to the semielliptical 
surface crack problem works in the following way: 
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1. The shear loading, t zx or T Zy > is applied to the crack surface 
as prescribed from the statement of the particular problem 
being solved. Solution 1 is used to compute the normal and 
shearing stresses at points on a plane in the infinite solid 
at the location of the front surface* 

2. The front surface of the body must be free of stresses, so the 

stresses produced in the previous step must be removed. This 

requires that opposing surface tractions be applied on the 

plane of an uncracked half-space. This is done by dividing 

the front surface into a number of rectangular areas and 

applying stresses to each area which are equal in magnitude 

but opposite in sign to the stresses computed at the center of 

each rectangle in Step 1. Solution 2 is then used to compute 

the shear stresses x and x at 96 points at the location of 

zx. zy 

the crack surface in the hal space. The effects of all 
surface rectangles are included by summing over the entire 
number of rectangles used on the surface. The result of this 
is a secondary residual stress which is now generated on the 
crack surface. 

Figure 2 shows a typical arrangement of surface 
rectangles used in the analysis. This pattern is the one 
which appears in the first quadrant and is repeated in each 
of the other three quadrants. 

3. The stress residuals which have now been produced on the crack 
surface must be removed in order that the boundary conditions 
on the crack surface can be met. This is done by applying 
opposing stresses on the crack surface. The polynomial 
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i 

i 

equations (1.3) for and x^ are least square fitted to I 

the values of the stress residuals at the 96 points on the 
crack surface shown in Fig. 3. 

I 

The change in the stress intensity factors AK ? and AK 

i s t 

caused by the application of these opposing stresses is 
calculated and added to the stress intensity factors caused by 

1 i 

the initial loading. Solution 1 is again used to compute the | 

j 

stresses on the front surface if more than one iteration is 

1 i 

I 

desired on that surface. Either one or two iterations of this | 

kind are done with the front surface, storing the stress I 

i 

residuals applied to the crack surface as well as the stress j 

residuals applied to the front surface. | 

4. When stress residuals are applied to the front surface, 
stresses are produced at the location of the back surface. 

These stresses must also be accounted for during the applica- 
tion of the alternating technique. The applied stresses on 
the front surface are used to calculate the stresses on the ■ 

t 

back surface, using Solution 2, and stresses of the opposite 
sign are then applied to the back surface. This, in turn, 
causes some additional stress residuals on the crack which 
must be removed, and the change in the stress intensity 
factors AK., and A is once again calculated. This change is 
then added to the running sum of each stress intensity factor 

t 

which was calculated from the initial loading plus the contri- 
butions of any subsequent iterations. 


5. The same steps which have been accomplished between and crack 
and the f-ont surface are now done between the crack and the 
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back surface, and a process similar to the one between the 
front and the back surface is now performed between the back 
and front surface. 

Steps one through five constitute one cycle of iteration. Several 
cycles are executed until the changes in 1^ and from each step are 

negligible as compared to the values caused by the initial loading. 

The final stress intensity factors are the sum of all the stress 
intensities due to the initially prescribed crack surface stresses plus 
the contributions due to all the stress residuals which arc applied at 
the crack surface. 

This procedure is extremely complicated and time-consuming when it 
is considered that a large number of rectangles are used on each 
surface, and 96 points are used on the crack face. For this reason, 
all such calculations were programmed and performed on the digital 


computer . 
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Chapter 2 

THE ELLIPTICAL CRACK WITH NONUNIFORM SHEAR LOADING 


1 .0 The Problem and its Boundary Conditions 

The problem, referred to earlier as ’’Solution 1", is illustrated 
in Fig. 4, and is a flat elliptical crack embedded in an infinite 
solid. The crack has a semiminor axis "b" located along the y-axis, 
and a semimajor axis ”a” located along the x-axis. The plane of the 
crack is considered to lie on the plane defined by z - 0. The surface 
of the crack is subjected to nonuniform shear stresses and the 

positive directions of which are shown in Fig. 4. 

Considerations of the conditions of skew symmetry of the displace- 
ments lead to the following boundary conditions: 


and 



where u and v are the displacements in the x and y-directions, 
respectively. It is assumed that the shear stresses applied to the 
surface of the crack are expressible in terms of a polynomial of x and 
y as 


zx 


EE 

m=0 n=0 


* m n 
A x y 
mn 7 





3 3 

EE 

m*0 n=0 


o m n 
B x y 
mn 7 


(2.3) 
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where m + n 3. 

All stresses are assumed to vanish at large distances from the crack. 

2.0 Stress Functions, Displacements and Stresses 

Kassir and Sih (22) formulated the problem of an elliptical crack 
subjected to uniform shear applied in any direction as a mixed boundary 
value problem. They were able to satisfy the three-dimensional equa- 
tions of elasticity through the appropriate selection of two harmonic 
stress functions, f and g. These functions satisfy Laplace's equations 


v 2 f = v 2 g = 0. (2.4) 

Navier's displacement equations of equilibrium are satisfied by the 
harmonic functions, f and g. Define 

F * |£ + |* and V 2 F = 0. (2.5) 

ax ay 

The rectangular components u, v, and w of the displacement vector 
are given in terms of the harmonic functions as: 

u = -2(1 - n) ft ♦ z (2.6a) 

v = -2(1 - n) |f ♦ * ff (2.6b) 

w = -(1 - 2n) F + z (2.6c) 


where u, v, and w are the displacement components in the x, y, and z 
directions, respectively, and n is Poisson's ratio. The strain- 
displacement equations may be applied to the displacements to yield the 
strains, and in turn the stress- strain relationships for linear 
elasticity may be used to obtain the stresses in terms of the stress 
functions f and g to be 
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-** = - 2(1 
2G 


, » 2 f , 9F S 2 F 

n) 5151 • 2n 3? * 1 77 

3x 
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-£ 


- 2(1 - n ) 


3y3z 


0 3F 3 2 F 

' 2n ^* v 


°zz 3^F 
2G = 2 3z 2 


~xr = _(1 - n) ^ (if + lf) + 


3 2 F 

3x3y 


xz 

2G 


2G 


3 2 f 

3F 

- 

3z 2 

n 3x 

3x3z 


3F 

n — + 

3 2 F 

3z 2 

n 3y 

s 3y3z 


(2.7a) 

(2.7b) 

(2.7c) 

(2 . 7d) 
(2.7e) 

(2 . 7f) 


where G is the shear modulus. 

It has been found convenient in solving boundary value problems of 
this nature to introduce the ellipsoidal coordinate system X, u, and v 
where these ellipsoidal coordinates are roots of the equation 

2 


1 - 


a +s b +s 


2 2 

f-- 0. 

2 s 


( 2 . 8 ) 


The ellipsoidal coordinates X, y and v are related to the cartesian 
coordinates x, y and z by the transformation equations: 


2 2 2 2 2 2 2 

a (a -b ) x * (a z +X) (a%y) (a +v) 

2 2 2 2 2 2 2 

b (b -a*) y^ = (b +X) 0> fy) 0> +v) 


(2.9a) 

(2.9b) 

(2.9c) 


2 ■> 2 

a b z = Xyv 

where the values of the coordinates are limited to the standard 
restrictions : 

2 2 

-a < v < -b <y<0<X<«° (2.10) 
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In these coordinates, on the plane z = 0, the interior of the crack is 
described by X = 0 while the region exterior to the crack is described 
by v = 0. 


Segedin (25) suggested the use of a series of functions which 
satisfy (2.4) and have properties which make it possible to satisfy 
boundary conditions of the type required by (2.1), (2.2), and (2.3). 
Shah and Kobayashi (6) used them to develop the solution for an 
embedded elliptical crack subjected to nonuniform normal loading. 
These functions are used here in a similar way to solve this problem 
and are expressed in the form 


f \ 3 3 ' a mn\3 m+n v (m+n+1) 



'EE 


g/ w.0 wO \ b mn / 9x " ^ 


where m + n < 3, 


-o> . f Jjh ds 


x 'Q(s) 


( 2 . 11 ) 


( 2 . 12 ) 


u»(s) 



(2.13) 


and 

Q(s) = s(a 2 +s) (b 2 +s). (2.14) 

Segedin (25) showed that is a harmonic potential function. Since 

V (n) is harmonic, the stress functions, f and g, are also harmonic and 
therefore satisfy the governing differential equations (2.4). 


3.0 Satisfying the Boundary Conditions 

The boundary conditions, as defined in (2.1), (2.2), and (2.3), 
must now be satisfied to complete the formal elastic solution. 

Equation (2.1) states the first boundary condition as being a * 0 on 
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the plane z = 0. It is obvious from (2.7c) that this boundary condition 
is identically satisfied at any point on the plane z = 0. 

Equation (2.2) states that the displacements in the x and y-direc- 
tion, u and v, are equal to zero outside the ellipse on the plane z = 0. 
The substitution of this equation into (2.6a) and (2.6b) yields: 


if = 

3z 3z 


0 


(2.15) 


Segedin (25) showed that the z -derivative of the function vanishes 

on the plane z = 0 in the region outside of the ellipse. From (2.11) 
it can be seen that if 3V^ m V<>z vanishes outside the ellipse for z = 0, 
then 3f/3z and 3g/3z must likewise vanish, thus satisfying the boundary 
conditions (2.2). 

The final task remaining in satisfying the boundary conditions is 

that of satisfying (2.3) for and x inside the crack on the plane 

z = 0. This was done by deriving the expressions for the stress 

components and t ^ on the plane z = 0 inside of the ellipse, using 

(2.7e) and (2.7f), for each term in the polynomial which is produced by 

(2.11). The result was a polynomial expression corresponding to each 

one of the ten independent a coefficients and the ten independent b 

mn mn 

coefficients. To satisfy the shear stress boundary conditions, 
coefficients of similar powers of x m y n were collected from all the 
calculations described above and equated to the corresponding coefficient 
from (2.3). This resulted in the following matrix equation: 


is--"** I 




m 


v 



! 
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i 



In the column matrices, the terms b _ and B are permuted in an 

mn mn r 

identical manner to the a and A terms, which are shown as the top 

mn mn r 

half of the column matrices. The square matrix is a 20 x 20 matrix 

and nonzero terms of this matrix are tabulated in Appendix I. 

The solution to a typical problem would require that the A and 

®mn coe ^fi c i en t s the right-hand column matrix be known or specified 

for the particular problem being solved. It is quite likely that the 

problem being solved would require only a portion of these coefficients 

to be nonzero. Tt is also quite likely that the symmetry or skew- 

symmetry of the stresses being used could preclude the necessity for 

' using all twenty of the coefficients. Equation (2.16) may be solved for 

t the a and b_ coefficients. These would be substituted into (2.11) 

mn mn 

and then into (2.5) and (2.7) to calculate the stresses, and into (2.6) 
to compute the displacements. 
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This completes the formal solution to the three-dimensional 
equations of elasticity. 

4.0 Stress Intensity Factor s 

The mode two and mode three stress intensity factors would be o r 
more importance to the designer than would the stresses or the displa:e- 
ments near the crack. These stress intensity factors, when defined as a 
function of position along the crack border, would be used with a 
suitable failure criterion to help accomplish the design of mechanical 
parts loaded with nonuniform shear stresses. 

4.1 Definition of K 2 and 

Mode two and mode three stress intensity factors have been computed 

for each a _ and b term in the stress functions of (2.11). This was 
mn mn ' 

done by computing the following limits for mode two and mode three 
stress intensities, respectively: 



Figure 5 illustrates the terms used to define K 2 and in (2.17). The 
direction normal to the crack front is defined as n, while the direction 
tangential to the crack front is defined as t. The angle <|> is the 
parametric angle of the ellipse, and the angle 6 is the angle between 
the x-axis and the outward normal vector. It should be noted that the 
angle 4 locates a point along the border of the ellipse which does not 
in general coincide with the intersection of the radial line with the 
ellipse. In these expressions, r is the distance from the crack tip 
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measured along an outward normal, and t and t refer to the stresses 

zn z t 

in the normal and tangential directions, respectively. 

4 . 2 Derivations of Stress Intensity Factors 

Mode two and mode three stress intensity factors, and K^, are 

calculated for every coefficient a and b . In the following para- 

mn 

graphs a calculation of the stress intevs:*- • factors corresponding to 
the coefficient a^ will be done to illu_' .-..Le the procedure. 

The first step is to compute the stress components and x on 
the plane z = 0 outside the crack. liquation (2.11) is substituted into 
(2.7e) and (2.7f) to calculate these stresses. The use of Leibnitz' 
rule for differentiiting integrals leads to: 


2G T zx 


z«0 


-8(1 - n) a 1Q x r 4 


* 8na in x r 5 
z*0 10 5 lz-0 

p=0 y*0 


(2.18) 


2G T zy 


z-0 


8na lO y P 6 


2*0 

M*0 


(2.19) 


where the I\ -terms are shorthand notation for: 


r, = 


on 

f ds z 3_x_ 

a S (a 2 .5)^OT " A(a 2 .X)/5W 51 


( 2 . 20 ) 


09 

■ \I 77 


r m , I ds x 3X_ 

5 i .J , 2 % 2 J~T ' 7 2 ,'t_j~7rr 3x 


X-' <.%.)‘-ffiO (aW»W 


( 2 . 21 ) 


0 

•/ 


ds 


2 ,.. t 2 


A' (a%s)(b‘+s)/Q(s) (a 2 *X)(b 2 *A)/Q(x) 


ii 

dx 


( 2 . 22 ) 


Byrd and Friedman's (26) book is used to calculate the above 
intcgrnls directly in terms of elliptic integrals and Jacobian elliptic 
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functions. The derivates of X may be obtained by differentiating 
(2.9a), (2.9b), and (2.9c). Equations (2.20), (2.21), and (2.22) 
become 

r 4 * ~i ((>-2k 2 )E(u : ) * Anu^nu, - V 2 u,] 

(2.23) 

2 2 

. 2z (b +X) 

X/Q(XT 0 4) (X-v) 

I'j * |(2*k 2 )Uj-2(l*k 2 )E(Uj)+k 2 snU|CnUjdnUjj 

(2.24) 

_ 2x 2 X(b 2 +X) 

(a : + x) 2 /QTxT (X-u) (w-v) 


6 aYV 


(l+k' 2 )E(u 1 )-2k ,2 u 1 -k 2 snu 1 cdu 1 


2x 2 X 


(a +X) rQ(X) (X-u) (X-v) 


(2.25) 


where E(u^) is the incomplete elliptic integral of the second kind. The 
modulus k and the complementary modulus k' are defined by 

k 2 « 1 - b 2 /a 2 ; k' - b/a (2.26) 

The functions sn, cn and dn are Jacobian elliptic functions and u^ 
represents the incomplete elliptic integral of the first kind. To find 
the limiting form of the stresses near the crack tip it is necessary to 
set u 3 0, and approach the crack tip. As r, the distance from the 
crack :ip, approaches zero the following limiting expressions are 
applicable: 



] 


5 


t 
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x -► 2abr 4> ; 


u = o ; 


v -♦ 


snuj -*■ 1 ; 


dnUj -*■ b/a ; 


cnuj -*■ 0 ; cdUj -*■ 0 (2.27) 


tnUj •+• a/X ; 


Uj •+ K(k) ; x -*■ a cos4> ; 


y -*• b sin<|» 


In these expressions 4> is given by 


* 2 . 2 , .2 2 . 
Q = a sin <t> + b cos <J> 


(2.28) 


and K(k) is the complete elliptic integral of the first kind. Upon 
substituting (2.27) into (2.23) through (2.25) and retaining only those 


terms containing (2r) the stresses are found to be 


^ T - = 16a io JlS? (' (1 ' n) cos *- n J7T cos "* 

3. K ° 


_i 

2 G 1 zx 


4> ' 


a ♦ 


) 


2G T zy * ' 16a 10 


nr / 

V£( 


n . . 2. 

si n<j> cos $ 


3/2. 
a $ 


•) 


(2.29) 


The stress components x and x . are now found from a stress 
r zn zt 

transformation to Le 


x * t cos 6 ♦ x sin 8 
zn zx zy 


x » -x sin 8 + x cos 8 
zt zx zy 


(2.30) 


The angle 8 is related to the parametric equations of an ellipse by 


. „ a sin* 
sm8 = r - 


„ b cos$ 
COS 8 * r— 2 * 

0* 


(2.31) 


Substitution of (2.31) into (2.30) yields the following equations for 


x and x „ . 
zn zt 
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T = T 

zn zx 


b cos<J> 


a sin<fr 


zy ** 


x - -t + T L£jsi 

zt zx ^ zy ^ 


( 2 . 52 ) 


Substituting (2.29) into (2.32) and the results of this into 
(2.17) yields: 


K. 


2 


K 


3 


-32G* 5 -J* .b. 2 X 

372 * ^ *10 C0S ♦ 

(ab) ' ® 10 

(2.33) 

32Gir* (1-v) A 

7T37r * *10 sin * cos * 

(2.34) 


(ab) 


Similar calculations have been done for every term a and b . 

ion mn 

The resulting equations are presented in Appendix II. 



Chapter 3 

THE HALF-SPACE PROBLEM 


1 .0 Discussion of the Method 

This solution, which was referred to earlier as "Solution 2", gives 
the stresses in a half-space subjected to uniform normal and shear 
tractions over a rectangular portion of the surface of the half-space. 
The general formulation for this problem was originally developed by 
Boussinesq and Cerruti (23) for the cases in which the rectangular area 
is subjected to surface tractions normal and parallel to the surface of 
the half-space. Love (24) worked the problem in detail for the case in 
which the rectangular area is subjected to uniform normal stress. Alavi 
(11) gave the solution for the case where the rectangular area is 
subjected to shear tractions. 

Love's solution (24) gave the six stress components caused by the 
application of a load normal to the plane. Two of these components, yx 
and yz in Love's notation are needed for Solution 2. However, no 
derivation for these two components existed for the cases in which the 
two shearing components of the traction vector were applied, so it was 
necessary to derive these equations. 

The equations for the x, y, and z displacements are presented by 
Love (23) as: 


j_r 

2m [ 


3F 

1 

3H 

X 

3t|»j 

i 

dip 

(3 11 

3z 

2 

3x 

2(X+u) 

3x 

2 Z 

3x. 


3G 

3z 
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2 

3H 

ay 

A 

* 2(X+pT 

3y 

-h 

3*1 

3yJ 
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where X and u are Lame's constants defined as 

i - En 

(i+n)(i-2r0 

w * shear modulus 

and 

E * modulus of elasticity 
n = Poisson's ratio 

The harmonic functions F, G, H, and are defined in Love 
(23), and Smith and Alavi (12), and will not be further defined here 
since this solution may be explained with no additional knowledge of 
them. 

Equations (3.1), (3.2), and (3.3) are written in terms of what 
shall hereafter be called Love’s coordinate system, which is different 
from the one used in Chapters 1 and 2. This coordinate system is used 
in Figs. 6 and 7 where it may be seen that in terms of Love's coordi- 
nates the desired shear stresses on the surface of the crack are x 

y* 

and x . Love's notation will also be used to specify which component 
yz 

of the surface traction vector is being considered. In this notation, 

X^, Y v> and represent the applied load at the center of each rectangle 
in the x, y, and z directions, respectively. The derivations involving 
these components will be considered separately, and the derail solution 
will be obtained by superimposing the results from each. 

2.0 Derivation of Stresses on the Crack Surface 

The strain-displacement equations can be applied to (3.1), (5.2), 
and (3.3) to provide equations for the strains. Then the stress-strain 
relations can be used to yield equations for the desired stress 
components, x yx and x . These steps are accomplished in the next three 
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sections for the application of X^, Y , and to a ty-pinl rectangular 
area in the first quadrant of the plane. 

2.1 Stress Derivations when Applying X v 

figure 6 shows the coordinate system for a t/pic^ r^'-tangle used 
in this half-space solution. A shearing traction vector X v is applied 
parallel to the x-axis on the surface of the rectangle. This is 
described by the following boundary conditions. 


T zx = * X V < a * < b ’ z * 0) 


o =0 =0 =t = t =0 (lx I > a or lyl > b) 

xx xy zz xy zy v ' 1 w 1 


(3.4) 


Calculations made for t using the steps mentioned in the previous 

yx 


section yield: 


where 


T yx = 2 * [ 


3 2 F A 9 *1 

3y3z + (A+y) 3x3y " Z 3x3y 


3 2 ij) 


(3.5) 


3 2 F 


3y3z ~ X \> loR Taj+fa-x)’ 


d 4 ^(a-x) 


b 2 -(a+x) 


c 3 -(a+x) 


(3.6) 


3 2 ^ 

3x3y 


.. f a-x a»x a+x a-x 1 . 

= v L a i( z+a P b 2 (z+b 2 ) + c 3 (z+c 3 ) 


7) 


***i . * 

2 W x 


\— 

” L a i 


-X + a + X a + x _ a *X 

b 2 °3 d 4 


] 


(3.8) 


The terms a^ b,, c 3 , and d 4 are defined as 

= [(a-x)“ ♦ (b-y ) 2 ♦ z 2 J 

r i -> '>1 ' s 

b 2 * (a+x)~ ♦ (b-y) “ + z I 


.^iLa * — 


c 3 1 [( a+x ) 2 * (b + y) fc ♦ z 2 ] 

d 4 = [(a-x) 2 + (b+y) 2 + z 2 ] 


Similar calculations made for the other shear component on the 


surface of the crack give 


■±[- 


a 2. ' 
i Ul 
9y9z 


(3.9) 


where 


3 2 4» _ x 

9y9z ~ v 


r_i_ j_ _i_ j_' 

L a l ' b 2 + C 3 " d 4. 


(3.10) 


2.2 Stress Derivations when Applying Y v 

A shearing traction vector Y y is applied parallel to the y-axis, 
with the following boundary conditions: 


T zy * - Y v < a * lyl < b * 2 = °) 


(3.11) 


a =o = o = t = t * 0 , (lx > a or y > b) 
xx yy zz xy zx 11 11 


The following equation for t is obtained: 


__ 

yx 2 tt 


9 2 G A 3 *1 _ 9 2 » 

3x3z (\+y) 3x9y 3x9y 


where the derivatives were found to be 


(3.12) 


3 2 G y , 
9791 " Y v l0 « 


d 4 -(b+:J b 2 +(b-y) 

a^b-y) Cj- (b+y) 


3 i * y . b -y + b -y + _b n 

9x9y v a^(z+aj) b 2 (z+b 2 ) Cj(z+Cj) 


Lh = z JJL ♦ y 

9x9y 9x3y v 


hzl _ biil . b+jr + b+ 




(3.13) 


(3.14) 


(3.15) 
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Calculations made for t give 

yz 6 


yz 2tt 


d 2 c. 

3z 2 


a 2 j» 

3y3z 


where 


JLjl = 

3y3z 




a-x a+x . 

,2 2 1 a, b. 

(b-y) + z \ 1 2 


b+y [ a-x a+ x 
(b+y) 2 +z 2 \ d 4 C 3 


(3.16) 


(3.17) 


G is a harmonic function, so 


3 2 G 

3z 2 


'2 2 ' 
3^ 3 G 

2 + 2 
3x^ 3y' \ 


(3.18) 


3 2 G y 
3x 2 " V 


tan 


-1 b-y _ -1 b+y ^ -1 z(b-y) ^ -1 z(b+y) 

— + tan — *- - tan — (-* tan — f4- + 


a-x 


a-x (a-x)a^ (a-x)d 4 


+ tan * ^ZjL + tan * - tan , ... 

a+x a+x (a+x)b 


1 Zib^_ _ -> 2(blXi 


(a+x)c, 


(3.19) 


3 2 C 

7 ” 


-1 a-x _ -1 a+x -1 z (a-x) 

tan r — + tan r tan /w ~ -r r 

b-y b-y (b-y)a 


1 


tan' 1 l l a .t * L + 

tan (b-yib 2 


-1 a-x ^ -1 a+x -1 z(a-x) -1 z(a+x) 

* tan tan 517 * tan (b+yTd^ ■ tan fb^ 


(3.20) 


2.3 Stress Derivations when Applying 

The normal traction vector z^ is applied parallel to the z-axis and 
perpendicular to the plane, providing the following boundary conditions: 

o * -z (|x|<a, |y|<b,z«0) 

ZZ \) 

(3.21) 

o *0 *T =T = t = 0 ( J x I > a or |y| > b) 

xx yy xy zx zy 
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Calculations made for t yield 

yx 7 


r 2 2 i 

u 8 H djp 

T yx ~ 2rr (A+u) 9x3y + Z 3x3y 


where 


9 2 h (z+a 1 )(z+c 3 ) 

Wy = Z v l0g (z+b 2 ) (z*d 4 ) 


(3.22) 


(3.23) 


l±, z i 

8x3y v a 


J_ + J_ _L 

1 b 2 C 3 d 4 


(3.24) 


The following equation is obtained for t 


yz 


yz 2w 


-z 


a 2. 

8 y> 

8y8z 


(3.25) 


where 


o 

8 " 4 » 

8y8z 


2 2 
(b-y)^ 


a-jc a+x 
a l + b 2 


2 2 

(b+y)%z z 


a-x 

d. 


♦ 

C V. 


(3.26) 


3.0 The Solution Form for Application to the Crack Problem 

When using the half-space solution, each quadrant of the plane is 
divided into a number of small rectangles, as shown in Fig. 2. The 
exact number of rectangles used can vary, depending upon the particular 
crack shape a/2c, and the amount of crack penetration into the plate, 
a/t. The effect of applying the normal and shear loading to each 
rectangle must be considered when the stress at any particular point is 
calculated. The number of calculations would tend to get out of hand if 
a large number of rectangles were used, so it is desirable to take 
advantage of the symmetry conditions of the stresses and thus reduce the 


i 
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number of total calculations necessary. Therefore, calculations are 

made using only the rectangles in one quadrant of the surface. 

The principle of superposition is used to determine the shear 

stresses on the face of the crack, x and x , by adding the effects of 

yx yz ' 6 

applying X^, Y , and Z^. In order to clarify the procedure, and shorten 

the length of the equations used, the following shorthand notation will 

be adopted. Let (3.5), (3.12), and (3.22) be called FIX, F2X, and F3X , 

respectively. The number 1, 2, or 3 identifies with the direction of 

loading, X^, Y , or Z^, respectively, while the last letter X identifies 

that the shear component in the x-direction, t , is being determined. 

yx 

In a like manner, let (3.9), (3.16), and (3.25) be identified as F1Z, 
F2Z, and F3Z, respectively, where the Z denotes that the stress 
component in the z-direction, x is being considered. 

Using supexposition, and factoring the applied load from each of 
the equations mentioned above, the x component may be written: 

t * X FIX ♦ Y F2X ♦ Z F3X (3.27) 

yx v v v 

Consider the case in which a uniform load is applied to the surface of 

the crack in the x-direction only. This would correspond to the loading 

being applied in the direction of x as shown in Fig. 7, and would 

yx 

result in the stress symmetiies on the surface of the plane as shown in 
the same figure. These symmetries may be written as 
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in X, 

Odd 
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Odd 

in Y 


(3.28) 
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They will be used in generalizing the solution for computer use. 
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A check of the symmetries of the equations derived for FIX, F2X, 
and F3X, shows them to be: 

FIX = Even in X, Odd in Y 

F2X = OJd in X, Even in Y (3.29) 

F3X = Odd in X, Odd in Y 

To illustrate how these symmetries are used, consider only the 
effect of the X^ surface traction on four rectangular areas, located 
symmetrically in each of the four quadrants, as shown in Fig. 7. Then 


F IX j = FIX (X - X, Y - Y, Z), 
F1X 2 = FIX (X ♦ X, Y - Y, Z), 
FIXj = FIX (X - X, Y + 7, Z), 
F1X 4 = FIX (X + X, Y + Y, Z), 


first quadrant 
second quadrant 
third quadrant 
fourth quadrant 


(3.30) 


where X, Y, and Z are the coordinates of the point at which the stresses 
are being calculated, while X and Y are the coordinates of the center of 
the rectangle being considered. Since the plane is the origin of the 
coordinate system, the Z coordinate is equal to zero. On the plane of 
the crack, Y = 0, and after using -Y and the symmetry conditions of 
(3.29), equations (3,30) become: 


FIXj = FIX (X - X, -Y, Z) 
F1X 2 * FIX (X + X, -Y, Z) 
F1X 3 ■ -FIX (X - X, -Y, Z) 
F1X 4 = -FIX (X ♦ X, -Y, Z) 


From (3.31), FIX * -FIXj and F1X 2 « -F1X 4 . 
Let 


t = X FIX, ♦ X FIX. 4 
y* Vj 1 v 2 2 


( FIX, ♦ X 
v 3 5 V 4 


FIX, 


(3.31) 


(3.32) 
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where the subscripts on X y indicate the quadrant number. Using the 
symmetries of (3.28), 


X = x = -X = -X 
V 1 V 2 V 3 V 4 


(3.33) 


and substituting. 


t = 2 X I FIX(X-X, -Y, Z) + FIX(X+X, -Y, Z) 
yx v i i j 


(3.34) 


If a similar procedure is followed when applying Y^ and Z^, and if the 

results from all rectangles are summed, the stress t generated at a 

yx 


point on the surface of the crack is 


t = ♦ > 2 X 

yX k-1 V k 


£ 2 X u | FlX(X-X k , -Y k , Z) + FIX (X+X 


-k' -V Z >1 * 


+ £ 2 Y y |F2X(X-X k , -Y k , Z) - F2X(X*X fc , -Y R . Z) j ♦ (3.3S) 


♦ £ 2 Z JF3X (X-X , -Y fc , Z) - F3X(X+X k , -Y k , Z) 
k-1 k 1 


where k is a particular rectangle number and N is the total number of 
rectangles in the first quadrant. 

In the same manner as above, the t component of the stress may 

y z 

be written as: 


t = X FlZ ♦ Y 1-2Z ♦ Z F3Z 
yz v v v 


(3.36) 


The equations for F1Z, F2Z, and F3Z provide the following symmetries: 


F1Z * Odd in X, Odd in Y 
F2Z » Even in X, Even in Y 


(3.37) 


F3Z * Even in X, Odd in Y 
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If the same procedure is followed to calculate as was followed to 

cnlculate t in (3.35), the results will be: 
yx 


t = ♦ T 2 X 

y* k=i \ 


FlZ(X-X k , -? k . Z) ♦ FlZ(X+X k , -Y k , Z) 


E 2 Y 

k=l k 


F2Z(X-X k , -Y k , Z) - F2Z(X + X k . -Y k , Z) 


(3.38) 


E 2 Z 

k=l v k 


F3Z(X-X . -Y , Z) - F3Z(X+X k , -Y k> Z) 


Chapter 4 

THE COMPUTER PROGRAM 

1 . 0 General Makeup of the Program 

The mathematical derivations discussed in Chapters 2 and 3 result 
in an appalling number of equations which must be solved to provide a 
solution to the problem. Additionally, the alternating technique with 
its iterations causes many of these equations to be solved repeatedly 
during the solution procedure used for a single crack geometry. The 
magnitude of the problem dictates that the digital computer be utilized 
to obtain the solutions. Therefore, the entire problem was programmed 
into the Fortran IV language to be used on the Colorado State University 
CDC 6400 computer. 

Subroutines were used to a great extent in the makeup of the 
program. Generally, a separate subroutine was written to accomplish 
each individual portion of the overall solution. The subroutines were 
then called in the proper sequence necessary to successfully achieve the 
solution to the problem. 

As each subroutine which performed a particular function was 
written, it was checked and debugged individually or in conjunction with 
some other required subroutines. This facilitated the checkout of the 
large program when the subroutines were added to it since they had 
already been checked to be operational on an individual basis. After 
an initial period of debugging, the proven portions of the large program 
were put on permanent file in the computer. Five individual cycles of 
a single file name were used so that any changes in a particular sub- 
routine would necessitate changes only to one cycle, which constituted 
only a small portion of the total. 
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The following steps in abbreviated form show the operating sequence 
of the computer program. 

1. Re°.d the input data for the problem being solved. This data 

includes input loading in terms of coefficients A and R 

r 6 mn ron 

crack and plate geometries, number of cycle iterations to be 
made, and the dimensional informat ior for the rectangular grids 
used on both the front and back surfaces. 

2. Generate the coordinates for the points on the crack surface 
for which stress calculations are to be made. 

3. Generate the initial shear stress loading on each of the points 

on the crack surface from the input coefficients A . and B . 

r mn mn 

4. Perform a least squares fit of the stress distribution 
calculated in the last step. This step is not necessary when 
the stress distribution being fitted is the initial loading 
determined by the coefficients prescribed in step one, but it 
is done for this first iteration to keep the program general 
enough to handle successive iterations, and also as a check on 
th~ accuracy of the least squares solution for a known loading 
case. The method of least squ • s is discussed in more detail 
in Section 2.2. 

5. Calculate the modulus and complementary modulus of the Jacobian 
elliptic functions and integrals, and compute the complete 
elliptic integrals of the first and second kinds. All four of 
these functions arc used in the 20 x 20 matrix of 
equation (2.16). 

6. Using the functions of the last step, determine the nonzero 
term;? of the matrix. Use a Gauss elimination method to 


I 
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solve the matrix equation, and calculate the coefficients 


a and b 
mn mn 


7. Calculate the mode two and mode three stress intensity factors, 
K 2 and Nondimensionalize these factors and call them 

and M 3 , the stress intensity magnification factors. 

8. Calculate the stresses at the center of each rectangle on the 
front surface. This requires the following steps for each 
rectangle: 

8.1 Calculate the ellipsoidal coordinates X, u and v of the 
rectangle center. 

8.2 Calculate the required elliptic integrals, Jacobian 
elliptic functions, and partial derivatives. Compute 
the stresses at the center of each rectangle. 

9. Calculate the stresses on the crack due to freeing the stresses 
on the front surface. 

10. Perform steps four through seven for the new loading on the 
crack surface. Repeat steps eight through ten for as many 
iterations as are desired on the front surface. 

11. Obtain the total stress applied to the front surface and use 
this total to calculate the resulting residual stresses at the 
center of each rectangle of the grid on the bacx surface, as 
in 8.1 and 8.2. 

12. Calculate the stresses on the crack due to freeing the stresses 
on the back surface. 

13. Repeat steps four through seven for the new loading on the 
crack surface. 


| m mm 
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14. Obtair the total stresses applied to the track in steps three 
and nine, and repeat steps four through 14 for the back 
surface. This completes one cycle of iteration, and should be 
repeated until the stress intensity factors are negligible 
compared to the ones obtained from the original applied load. 

2.0 Special Subroutines Used 

Some of the subroutines used in Section 1.0 should be discussed to 
point out details of their operation. 


2.1 C rack Coordinate Generator 

A subroutine called COORD is used to generate the 96 points on the 
crack surface for which stress calculations are to be made. Considera- 
tions of symmetry enable calculations to be made on only the half of the 
crack for which x is positive, as shown in Fig. 3. 

The coordinate generator uses five similar ellipses, each with a 
semiminor and semimajor axis which is a particular fraction of the 
length of the axes of the outer ellipse. The ratios used in this case 
were 0.2, 0.4, 0.6, 0.8, and 1.0. Coordinates were generated by a 
simultaneous solution of the two equations 


and 


2 2 
x y _ , 

a b 


tan $ * y/ x 


(4.1) 


(4.2) 


where the angle $ was increased in 10 degree increments for each of the 
five similar ellipses. Values were calculated only in the first 
quadrant (positive x and y) and then reflected with the proper sign into 
the other quadrant. 
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2 . 2 Least Squares Curve fitting 

Tho loading on the crack surface for which all calculations were 

made in this w 'rk was one in which only a uniform x-componcnt t jx was 

applied as shown in fig. 4. A consideration of the problem shows that 

the t and t components applied back on the crack surface through the 
z x zy 

iteration process have symmetries which state that t zx is even in x and 

t is odd in x. Because of these symmetries, the least squares curve 
zy 

fitting process is slightly different for each of these two components. 
If t zx is even in x, the first of equations (2.3) becomes: 


T zx = A oo + A oi y * A ™ x + A n9 y ♦ A 9i x y + A r,x y 


20 


02 


*21 


03 


(4.3) 


The x zx stress components calculated at the 96 points mentioned in the 

last section are used in the least squares subroutine, LEAST, which 

calculates the fit for this component. 

The component t is odd in x, so the last of equations (2.3) may 
zy 

be written: 

V ■ 8 00 * B 10 x * B 11 xy * B 02 y2 * B 12 xy2 * B 30 x3 (44) 

After the t ^ components were calculated on the surface of the crack the 

method of least squares was used on the stress distribution to calculate 

the B coefficients in a subroutine called LEAST2. 
mn 

The least squares criterion used is a general one in which 

NT 1 2 

E - £ [P(x i ,y i ) - f(x..y.)J (4.5) 


where 

P( x i,yi) = tabular pressure values 

f(x it y^) * family of functions representing the pressure, which 
in this case would be (4.3) or (4.4). 
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The best fit will be obtained if E is a minimum. This is obtained by 
setting 

* 0 < 4 - 6 > 
mn 

where represents the coefficients in (4.3) or (4.4). 

2.3 Calculation of Elliptic Integrals 

The complete elliptic integrals of the first and second kinds are 
required in the elements of the matrix in equation (2.16), so they 
must be evaluated fcr their use there. A polynomial approximation, as 
presented in Abramowitz and Stegun (27) , is used to calculate the values 
of these integrals. The absolute value oi the error associated with 
these polynomials is no greater than 2.0 x 10 , so it is sufficiently 

accurate for the intended purpose. 

Incomplete elliptic integrals of the first and second kinds are 
required for the stress calculations of Solution 1. These integrals 
are evaluated using the process of the arithmetic-geometric mean, which 
is also presented in (27). 

2.4 Ellipsoidal Coordinates 

Ellipsoidal coordinates X, y and v were obtained as the roots of 
the cubic equation (2.8). The Newton-Raphson method was used to find 
the first root, and the quadratic equation solution was then used to 
« find the other two roots of the reduced polynomial. The three roots 

thus determined were placed in the order of ascending magnitude and 
each was matched with its proper ellipsoidal coordinate, as shown in 
(2.10). A final check was mode to assure that the coordinates did 
indeed have the proper relative magnitudes of (2.10), and could satisfy 
the transformation equations (2.9). 



3 . 0 Rectangular Surface Grids 

A typical rectangular surface grid used to calculate and remove 
stresses from the front and back surfaces is shown in Fig. 2. Most of 
the cases run utilized a grid with 62 rectangles on the front surface 
and a grid with 32 rectangles on the back surface. However, more 
rectangles were used in cases involving thin cracks or cracks which 
penetrated the plate for a major part of the thickness. The total 
computer time used is dramatically affected by the number of rectangles 
on each surface, so the number was kept to a minimum consistent with the 
desired accuracy and convergence required for each particular crack 
geometry. 







Chapter 5 


DISCUSSION OF RESULTS FOR PART I 

The results included in this section were obtained from computer 

calculations using the equations developed in the previous chapters. 

All results shown are for the loading case in which a uniform shear load 

was applied to the surface of the crack, paralljl to the x-axis, as 

shown in Fig. 4. The initial loading did not include the application of 

the shear component in the y-direction, t . This particular initial 

loading reduces the first of equations (2.3) to the case where all A 

inn 

coefficients except A Q0 are equal to zero, while all coefficients of the 
second equation of (2.3) are equal to zero. Solutions 1 and 2 of the 
previous chapters were used with the alternating method to determine the 
mode two and mode three stress intensity factors for a semielliptical 
surface crack in the surface of a finite thickness solid. 

1.0 Results and Discussion 

Mode two and mode three stress intensity factors were calculated 
for several semielliptical surface crack shapes subjected to uniform 
shear loading applied only in the x -direction. The various crack 
shapes considered had a/2c ratios of 0.05, 0.1, 0.2, 0.3, and 0.4, and 
had depth ratios a/t varying from 0.2 to 0.9. A Poisson's ratio of 0.25 
was used for all calculations. 

Figures 9 through 13 present the mode two and mode three stress 
intensity factors as a function of the parametric angle $ which is 
measured from the semiminor axis of the ellipse, as shown in Fig. 8. 

The stress intensity factors K 2 and have been nondimensionalized and 
plotted as M 2 and M^, respectively, where 
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M 2 = £ U 

and 

K 3 

M = ^ 

3 t NQ)^ 


( 5 . 1 ) 


(5.2) 


The stress t is the initially applied loading, "a" is the semiminor 
zx 

u 

axis, and Q = E(k), the complete elliptic integral of the second kind. 

The process of iterating between the crack and the front and back 
surfaces of the plate causes the and curves to assume new shapes 
and values as compared to those obtained for the embedded crack. 

Several effects interact to influence the variation of the shape of 
these curves, and an inspection of the step-by-step output from the 
computer analysis reveals the following trends: 

1. The crack shape effect: This effect is a result of the 

elliptical crack solution. It causes to increase and 
to decrease as the angle 4> varies from 0 to 90 degrees. The 
range through K 2 increases becomes smaller as the crack shape 
goes from a/2c * 0.40 to 0.05, while the range through which 
decreases becomes larger with the same variation in crack 
shape . 

2. The front surface effect: This effect increases the value of 

Kj as the angle $ varies from 0 to 90 degrees, and is larger 
for cracks having large a/2c ratios than it is for the slender 
cracks having small a/2c ratios. It is this effect which 
produces a small value for at $ * 90° on cracks having 
a/2c = 0.20 and larger. 
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3. The back surface effect: This increases the value of by an 

amount, which is maximum at $ = 0° . As $ is allowed to range 
from 0 to 90 degrees the increase in diminishes and 
essentially disappears at about 60 or 70 degrees. This back 
surface effect also causes to increase from zero at 4> » 0° 
to a maximum at approximately 40 or 50 degrees, and drop back 
to a minimum value at d> = 90°. The magnitude of change in 

is substantially less than that in Kj, especially for a/t 
ratios greater than 0.5, where the magnitude of is 
substantial . 

4. Effect of the interaction between the front and back surfaces: 
This effect causes both K 2 and to increase. The increase 
of '<2 becomes greater as <J> varies from 0 to 90 degrees, while 
the increase of becomes smaller as $ ranges from 0 to 

90 degrees. 

It should be noted that the curves in the figures are all 
concave downward, while the K 2 curves are concave downward for the 
larger values of a/2c and concave upward for a/2c values less than 0.2. 
The results of this study, and those of Kassir and Sih (22) show that 
for an ellipse embedded in an infinite solid with only the t jx component 
applied, the curvature of the plot changes sign at an a/2c value of 
approximately 0.2. The behavior is a result of the crack shape effect 
for the particular loading used. 

Figures 9 through 13 show that the values for at 4 * 0 degrees 
increase with larger depth ratios a/t. This is due primarily to the 
back surface effect which exerts its maximum influence at ♦ » 0" . This 
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effect becomes increasingly important as the crack nears the back 
surface with large a/'t ratios. 

The figures will show that cracks with higher a/t ratios exhibit 
greater values than do those with lower a/t values. This is true for 
every value of 4> except 0 degrees where for all cracks is zero. 

These increases with higher a/t values can be attributed to the back 
surface effect, and the effects of the interaction between the front and 
back surfaces. 

All depth ratios a/t provide a value of zero for at 4> * 0°. 

According to equation (2.17) can exist only if the stress component 

T zn exists at this point. At the angle $ » 0° the tip of the crack is 

perpendicular to the y-axis, so the T zn component in this case is t 

which is an odd function of x. Since x * 0 at $ * 0°, then t must 

zy 

also be zero at this point, and K 2 must be the same. This justifies the 
fact that K 2 at = 0° is always zero, and also provides a check on the 
computer results for that point. 

Figures 9 through 13 show that the stress intensity factors 
approach a finite value at the angle * 90°. It is felt, however, that 
the curves may not represent the true behavior of the elastic solution 
where the crack tip intersects the surface of the plate because it is 
expected that there is a change in the nature of the crack tip singular- 
ity at that point. Hartranft and Sih (28) used a refined numerical 
analysis to study the problem of a semicircular surface crack subjected 
to normal loading. Their analysis utilized a very large number of 
surface rectangles and terms in the mathematical expressions used to 
calculate the stresses, and considered the singular and nonsingular 
portions of the surface stresses separately. They noted that the stress 
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intensity factor for the problem which they solved would tend to 
approach zero at the surface, though they were unable to extend the 
calculation to the point of intersection. Smith (29) concluded that the 
results obtained by Smith (8) and Thresher (14) agree well with those of 
Hartranft and Sih over most of the crack front although the smaller 
number of rectangles and mathematical terms used by Smith and Thresher 
miss the effect of the decrease in the stress intensity factor near 
« - 99°. 

It is expected that a comparable uncertainty exists in the behavior 
at the surface of the plate for the shear problem. The analytical 
approach used in this study, however, should be accurate along most of 
the crack front. The value of the stress intensity factors shown in the 
curves is representative of the average near <J> * 90°. 

The cost of the computer time necessary to perform a refined 
analysis such as the one by Hartranft and Sih becomes prohibitive. A 
serious study becomes necessary in which the advantages of greater 
accuracy are considered and weighed against the increased computer costs 
incurred. 

Figure 14 presents the values of maximum plotted against depth 
ratios a/t for cracks with a/2c values of 0.05, 0.1, 0.2, 0.3, and 0.4. 

It should be noted that the curve for a particular crack shape increases 
gently and smoothly as the depth ratio increases. There is a sizeable 
difference in the values for slender cracks with a/2c = 0.05 as 
compared with thicker cracks with a/2c values of 0.4. 

The maximum value of plotted in Fig. 14 is determined at 4 * 90° 
from Figs. 9 through 13. Considerations advanced in preceding para- 
graphs indicate that the maximum true value for may not exist at 
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(J) = 90° , Init this value is representative of in an average sense 
near <J> = 90°. 

Values of Kj for cracks with a/2c values ranging from 0.0S to 0.4 
are plotted against depth ratios a/t in Fig. 15. It should be noted 
that these curves for the corresponding range of a/2c values are 
clustered closer together than are those for K_. These K. curves also 
rise quite rapidly as the a/t ratio is increased hoyonci "pp'.'oximately 
0.6. This trend was not noticed with the curves, and indicates that 
cracks of this type are quite mode three sensitive to crack depth 
ratios . 

The results reveal that there is only a slight change in the 
maximum Kj and values for semielliptical cracks having low a/t ratios 
of 0.2, r Mnpared to the maximum and values for embedded 
ellipses of the same shape. The maximum value for a crack having an 
a/2c ratio of 0.2 increased by 15.1 and 31.2 percent as a/t was 
increased to 0.8 and 0.9, respectively, while the maximum K- value 
increased by 6.0 and 7.1 percent, respectively, under the same condi- 
tions. A crack having an a/ 2c of 0.3 showed an increase in Kj of 11.2 
and 22.6 percent and an increase in Kj of 5.3 and 6.7 percent as a/t was 
increased to 0.8 and 0.9, respectively. Similar increases were noted 
in both K 7 and for other crack geometries. 

2.0 Accuracy and Convergence 

Rectangular grids similar to the one illustrated in Fig. 2 were 
used on the front and back surfaces in the alternating technique 
employed f^.r this problem. It had previously been determined from other 
problems that convergence difficulties could exist if the rectangle 
sizes were too large for the problem geometry being considered. It was 
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discovered in the shear problem that convergence was not as difficult 
to attain as in the normal loading problem. The shear problem for a 
particular geometry could be worked using coarser grids and generally 
fewer cycles tha?> what were required for the normal problem. Usually a 
front grid with 62 rectangles, a rear grid with 32 rectangles, and two 
cycles of iteration were adequate to achieve excellent convergence for 
most problems. However, it was necessary to use four cycles and finer 
grids for problems having high a/t ratios, especially for those having 

I 

low a/2c values. 

Figure 16 shows the stress intensity factors after each cycle of 
iteration for a typical run in which four cycles were made on a crack 

i 

geometry of a/2c =0.2 and a/t = 0.9. It should be noted that the value j 

of the stress intensity factors increases fo” each cycle, and that the I 

i 

contribution of each successive cycle becomes smaller until the last 
cycle produces essentially an insignificant change. The differences 
between cycles three and four for are too small to be plotted on the 
scale used in the figure. 

The following tabulation shows the percentage of change of the 
maximum values of K 2 and K.. accomplished by each cycle of iteration. 

Cycle Percent Change in Percent Change in 

1 4.74 24.30 

2 1.42 4.09 

3 0.57 1.07 

4 0.22 0.35 

Obviously some truncation error exists in the final results due to 
stopping the iteration process ot a finite number of cycles. However, 
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enough cycles of iteration were performed in each case to insure that 
the change in the stress intensity factors which would be obtained by 
further iterations was of an insignificant magnitude and did not warrant 
the expenditure of additional computer time. 

Accuracy of the solution was checked by making two computer runs 
each on several cracks, using different surface grids for each of the 
runs on a specific crack geometry. The results of the runs for a 
particular crack configuration were then compared to determine the dif- 
ference which existed between them. Cracks considered in this manner 
had a/2c values of 0.05 and 0.1, and a/t values of 0.8 and 0.9. Solu- 
tions obtained from grids with 62 front and 32 back rectangles were com- 
pared with those ootained from grids with 55 front and back rectangles. 
Changing from one set of grids to the other caused the maximum values 
of to differ by 0.1 to 2.0 percent, and the maximum values of 
to differ by 0.6 to 0.9 percent. 

The results discussed in this chapter are obtained by numerical 
methods and thus include the step-by-step accumulation of numerical 
inaccuracies inherent in methods of this type. Although these results 
cannot be considered as an exact solution to the problem, they represent 
the best solutions available. 

The program was compiled under the Fortran Extended option which 
would give the minimum running time. It was then cataloged on per: t 

file to insure that no compiler time was necessary for the data produc- 
tion runs. Typical running time used by the CDC 6400 computer for 62 
rectangles on the front and 32 rectangles on the rear was approximately 
510 seconds for two cycles of iteration, and approximately 1015 seconds 
for four cycles. The cases in which 55 rectangles on each surface were 
used with four cycles required 1300 seconds. 
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PART 11 --NORMAL LOADING OF A SEMIELLIPTICAL SURFACE CRACK 

Chapter 6 

THE NORMAL LOADING PROBLEM 


1 . 0 Statement of the Semielliptical Surface Crack Problem 

Consider a flat semielliptical surface crack located in the ir^nt 
surface of a finite thickness slab. The geometrical configuration is 
the same as that shown in Fig. 1, but there is no shear loading applied 
to the surface of the crack for this problem. Instead, the surface of 
the crack is subjected to nonuniform normal stresses o 

z z 

The boundary conditions to be considered in the solution of this 
problem are as follows: 

T zx B T zy = 0 * (z * 0) C6 '^ 


and 


( 2 2 \ 
w “ 0 , l^T +y T >1 ' z ” 0 / 


(6.2) 


where w is the displacement in the z-direction. 

It is assumed that the normal stress applied to the surface of the 
crack is expressible as 


L L x y 


zz mn 

m«0 n=0 




(6.3) 


where m ♦ n <_ 3 because of the magnitude of work required in the ca’cu- 
lations if the sum were allowed to be larger than three. The front arid 
back surfaces of the body are considered to be free of stresses. 

2 . 0 The Method of Solution 

The method of solution which will be used for this problem is 
essentially identical to the one used for the shear problem which has 
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been discussed in previous chapters. The alternating technique is once 
again used with Solutions 1 and 2, defined as follows: 

Solution 1. The determination of the stresses near an elliptical 
crack at any point in an infinite solid through which 
it is desired to pass a plane. The surface of the 
crack is subjected to nonuniform normal loading. 

Solution 2. The determination of the stresses within a half-space 
subjected to uniform normal and shear stress over a 
rectangular portion of the surface. 

Solutions one and two are discussed briefly in the following two 
sections . 

2 . 1 Elastic Solution of a Crack in an Infinite Solid 

Solution 1 represents the elastic solution of an elliptical crack 
in an infinite solid subjected to nonuniform normal stresses. The 
particular form of this solution which was used for this portion of the 
problem was derived by Shah and Kobayashi (6). Thi method was used by 
them in an alternating method solution to the problem of an elliptical 
crack approaching the surface of an semi-infinite solid subjected to 
uniform tension perpendicular to the plane of the crack (17). Details 
and the derivations necessary to use the solution may be found in 
(0, 10. and 17), and will not be repeated here. 

2 . 2 The Half-Space Solution 

This solution has previously been identified as "Solution 2 ", and 
calculates the stresses in a half-space subjected to uniform normal and 
shear tractions applied to a rectangular portion of the surface of the 
half-space. The solution method is similar to that described in 
Chapter 5, and uses the equations descr ed in (11 and 12). 
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This portion of Solution 2 calculates only the o component of the 

stress on the crack surface. The symmetry conditions of the problem 

prohibit the presence of any residual shear stresses, t 7 or t , on the 

z y 

crack surface. Therefore, no attempt need be made at calculating then 
since they vanish when all effects from all rectangles are summed. 

2.3 The Alternating Method 

The alternating method used here is the same as the one used by 
Thresher and Smith (14, IS), and the method described in Chapter 1 which 
was used for the shear problem. It differs from the one used for the 
shear problem in that it calculates only the normal stresses on the 
crack surface and the resulting mode one stress intensity factor, , 
while the method described in Chapter 1 calculates both shear stresses 
on the crack surface and the mode two and mode three stress intensity 
factors, 1<2 and k,. 

As is the case with the shear problem, the alternating method 
iterations proceed until the change in the mode one stress intensity 
factor for any one cycle becomes negligible when compared to its value 
obtained for the initial loading. 

5.0 The Computer Program 

A computer program was written for the CDC 6400 digital computer 
which would perform all the calculations for Solutions 1 and 2. The 
general method of programming was the same as that discussed in 
Chapter 4 for the shear problem. In fact, the program for the normal 
loading problem was written before the one for the shear loading 
problem, and served as the model for it in addition to providing some of 
the subroutines which are all or partially common to both programs. 
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The sequence of operations performed and the manner in which the sub- 
routines are used is essentially the same ?r both programs. 

One of the primary differences between the two programs is the much 
smaller number of equations which must be solved for the normal loading 
problem. For example, the 20 x 20 matrix of equation (2.16) is only 
a 10 x 10 matrix for the normal load problem. Since only one component 
of stress is being applied to the crack instead of two, which is the 
case with the shear problem, there arc only approximately one -half as 
many terms necessary to calculate the stresses which must be freed from 
the surfaces of the body. 

A study of the symmetries of the normal loading problem reveals 

that the normal stress on the crack surface a is even in x. This 

zz 

allows (6.3) to be written as 
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(6.4) 


The method of least squares curve fitting is performed on this equation 
in the same manner as it was performed on the equations for the shear 
loading problem. After the curve fitting scheme is employed, a matrix 
equation similar to (2.16), but having a 10 x 10 K.j matrix, is solved 
for the ten independent a coefficients. These in turn are used to 
calculate the mode one stress intensity factor . 

A typical rectangular grid of the type used on the front and back 
surfaces is shown in Fig. I 7 , and will be discussed in mere detail in 


the next chapter. 



Chapter 7 


DISCUSSION OF RESULTS FOP. PART II 

All computer calculations were run for the case in which a uniform 

normal load was applied to the surface of the crack. This corresponds 

to the situation where a uniform censile load is applied to the body in 

a direction perpendicular to the plane of the crack. Equation (6.4) is 

then reduced to the case where all the coefficients A in the first 

mn 

iteration step are equal to zero except for A^. Solutions 1 and 2 were 
used with the alternating method to determine the mode one stress 
intensity factors for a semielliptical surface crack in the surface of 
a finite thickness solid. 

1.0 Results and Discussion 

Four different crack shapes were considered having a/2c ratios of 
0.05, 0.10, 0.20, and 0.30, and having depth ratios a/t ran .iig from 
0.20 to 0.95. All computations were made using a Poisson's ratio of 
0.25. 

The stress intensity factors plotted as a function of the paramet- 
ric angle $ of the ellipse are presented in a nondimensionalized form 
in Figs. IS through 21. The angle 4> is measured from the semiminor axis 
"a" at the portion of the crack which is farthest into the material, as 
shown in Fig. 8. Three effects interact to influence the variation of 
the shape of the stress intensity curves, and are: 

1. The effect of the crack shape: iTiis effect is a property of 

the elliptical crack solution and causes to decrease as the 
angle $ is varied from 0 to 90 degrees. The i.v ' of is 
constant for a circular crack with a/2c ratio of 0.5, out 
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decreases very drastically as <|> is increased for a slender 
crack with a/2c = 0.05. Thus the crack shape effect becomes 
more significant as the crack shape goes from a/2c - 0.5 to 
a/2c = 0.05. 

2. The effect of the front surface: Previous studies (9) and 

inspection of the present computer solutions show that this 
effect causes to increase as 4> varies from 0 to 90 degrees. 
The amount of increase from this effect is of smaller magnitude 
than the crack shape effect. 

3. The effect of the back surface: This effect can be considered 

to include the direct effect of the presence of the back 
surface as well as the effect of the interaction between the 
front and back surfaces. The back surface effect alone tends 
to cause a sharp increase in Kj near $ = 0°, while the 
interaction effect between the two surfaces tends to increase 
Kj near <J> = 90°, but by an amount substantially less than that 
caused by the back surface effect. The result of the two 
effects combined is to cause Kj to increase significantly near 
<j> = 0° and to increase by a smaller amount near 4> = 90°. 

The crack shape effect causes a general decrease in Kj with 
increasing $ , as can be noted in Figs. 18 through 21. The back surface 

effect causes an increase in Kj as the thickness ratio a/t is increased. 

In Figs. 20 and 21, for a/2c of 0.05 and 0.10, an increase in occurs 
before the expected decrease begins, which indicates that the maximum 
value of Kj does not necessarily occur at <p * 0°. The reason for this 

behavior is that for these slender cracks the front surface effect tends 



to increase Kj taster than the shape effect tends to decrease it in the 
range from 4> = 0° to approximately 4> = 40°. 

A cursory inspection of the curves of Figs. 18 through 21 reveals 
at least two general trends. It is apparent that for any given crack 
depth ratio a/t, the maximum Kj value is greater for slender cracks with 
small a/2c ratios than for the rounder cracks with larger a/2c ratios. 
Also, it is apparent that for any given a/2c ratio, the larger values 
occur for the deeper cracks, or those with higher values of a/t. 

The computations which were performed to obtain the results plotted 
in Figs. 18 through 21 used rectangular grids on the front and back 
surfaces similar to the one shown in Fig. 17. It was found that 
convergence difficulties existed for cases in which the plate thickness 
was approximately the same as the dimensions of the smallest rectangle. 
It was therefore necessary to revise the grid so that there were 
approximately twice as many rectangles near the origin of the grid. The 
revised grid has 79 rectangles, while the one shown in Fig. 17 has 55 
rectangles. The grid with 79 rectangles was used whenever the crack 
depth "a" approached the size of the smallest rectangle. Grids of 
62 and 32 rectangles were also used where the gee 'ry of the problem 
being solved would permit the use of these configurations. The use of 
the proper grid insured that good convergence was attained, as 
evidenced by the fact that the contribution of the last cycle was 
generally less than one percent of the total stress intensity factor. 

The stress intensity factor Kj , as shown in Figs. 18 through 21, 
smoothly approaches a definite value at <t = 90°. However, it is 
expected that there is a change in the nature of the crack tip singular- 
ity at the point where the tip of the crack intersects the surface of 
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the plate. This would mean that the curves do not represent the true 
behavior of the elastic solution there. Hartranft and Sih (28) have 
concluded that the mode one stress intensity factor tends to zero at the 
surface for a semicircular crack. They handled separately the singular 
and nonsingular portions of the stress on the surface of the plate. In 
addition, they used a very large number of rectangles and a large 
number of terms in the series expansions used to calculate the stresses. 
Smith (29) concluded that for the semicircular crack the level of 
approximation used here gives results within one or two percent over 
most of the crack border as compared with the results obtained from the 
more refined analysis of (28), although there can be a sizeable 
difference near the intersection of the crack tip with the front 
surface. 

It is therefore reasonable to expect a similar behavior for this 
analysis which irvolves a semielliptical crack rather than a semi- 
circular crack as considered in (29). The stress intensity factor 
at <)> = 90° has little practical importance, however, since the maximum 
occurs somewhere between the angle of 4> = 0° and 40°, as shown in 
Figs. 18 through 21 for semielliptical surface flaws. 

Figure 22 presents the values of maximum from this analysis for 
a/2c values of 0.05, 0.1, 0.2, and 0.5. The curves for a/2c values of 
0.4 and 0.5 were estimated by Shah and Kobayashi (17) from their solu- 
tion for the elliptical crack near the surface of a half-space. The 
results for the two-dimensional edge crack (50) which represents the 
limiting solution are also shown in the figure. It has been argued that 
cracks having an a/2c ratio less than 0.1 closely approximate the 
behavior of a two-dimensional edge crack, but these results indicate 
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that even a crack having a length 20 times its depth is somewhat less 
severe than a two-dimensional edge crack. 

2 . 0 Comparisons with Other Work 

A comparison of the results of this stud> with those due to Shah 
and Kobayashi (17) is made in Fig. 23. Shah and Kobayashi obtained 
their results by assuming that the back surface effects arc the same as 
the effect of the surface for an elliptical crack approaching the 
surface in a semi-infinite solid. They estimated the effect of the 
front surface, but chose to neglect the effect of the interaction 
between the front and back surfaces. 

The surface interaction effect causes only slight changes in the 
results for values of a/2c « 0.3 and greater. Though not shown for 
a/ 2c - 0.3, the results of this study compare well with those of Shah 
and Kobayashi for a/2c = 0.3, and they approach closely for a/2c = 0.2. 
However, the curves from this study become relatively higher as a/2c 
is decreased to 0.1 and 0.05, indicating the growing importance of the 
interaction between the front and back surfaces for slender cracks. 

The results of Fig. 22 were used to compute the fracture toughness 
values for a number of epoxy plates on which experiments were conducted 
by Smith (29). The fracture toughness was calculated as 

K IC = MjOfWQ) 5 * • (7.1) 

The nominal stress in the sample plate at the of fracture is 

indicated by a , and ”a" is the depth of the crack. The flaw shape 
parameter Q is defined as Q 2 * E(k), the complete elliptic integral of 
the second kind. The term Mj is the nondimensionalized form of Kj and 
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is frequently called the magnifi cion factor where the subscript 
indicates the mode one loading case. 

Figure 24 presents a plot of fracture toughness values which were 
calculated from the experiments of (29) using the curves cf Fig. 22. 

The data was used in equation (7.1) and was taken from a number of 
different epoxy plates having different average fracture toughness 
values. It was normalized by dividing the value for each sample by 
an average value for the epoxy plate from which the sample was 
taken. The average fracture toughness values which were used to normal- 
ize the data were calculated using only samples which had a/t values of 
0.5 or less. As shown in the figure, the a/ 2c values ranged from about 
0.13 to about 0.42. The normalized fracture toughness has a scatter 
band of approximately ±10 percent, and is relatively constant with 

variation in a/t, although it does drop slightly as a/t is increased. 

% 

The constant trend of the normalized values indicates that the curves 
are valid and possess some practical utility for design and analysis 
situations in which linear fracture mechanics may be applied. 

3.0 Advantages and Limitations of the Method 

The alternating method, as used in this study, is capable of 
producing a reasonably accurate solution to the types of problems 
discussed here. The method is straightforward and requires no approxi- 
mations from two - 'tnensional theories to effect the desired solution. 

F.xtremely good convergence can be realized but it is not always 
easily attained. As has been mentioned previously, convergence is quite 
dependent upon the number and size of the rectangles used in the grid. 
The desire for extreme accuracy would dictate the maximum number of 
rectangles possible. However, the total computer time used is a 
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function of the product of the number of front and back surface 
rectangles used, so a compromise is necessary. 

The total computer time required for some typical configurations 
for the shear loading problem has been discussed in Chapter 5. The 
normal loading problem was usually run with a total of four cycles of 
iteration for each crack geometry. This number of cycles required 
960 seconds of total computer time for 62 front and 32 rear rectangles, 
1220 seconds for 55 rectangles on each surface, and 2000 seconds for 
79 rectangles on each surface. The entire program was cataloged on 
permanent file within the computer, so none of the time mentioned above 
was lost to the CDC 6400 compiler. 

Accuracy of the infinite solid solution is dependent upon the 
preciseness with which the method of least squares is able to match the 
calculated stresses with the proper polynomial expression. A polynomial 
with more terms than the number allowed in a third- degree polynomial 
would undoubtedly increase the accuracy of the curve fit, but more 
terms were impossible to tolerate for reasons already mentioned. 

One of the limitations to the alternating method as it is presently 
used is its restriction to geometric configurations in which the front 
and back surfaces are parallel. The method is also limited to finite 
thickness bodies whose width is great enough so that no effects are 
realized from the presence of any side surfaces. This would pose a 
problem if solutions were desired for some geometries where the width 
effect could no longer be ignored. 
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Chapter 8 
CONCLUSION 


1.0 Discussion 

A semielliptical surface crack located in the surface of a finite 
thickness solid and subjected to two different conditions of loading is 
considered. The first condition of loading is one in which shear 
stresses are applied to the surface of the crack parallel to the semi- 
major axis. The second loading case is the one in which normal loading 
is applied to the surface of the crack in a direction perpendicular 
to the plane of the crack. The alternating method was used in each of 
the two loading cases to solve the problem. 

As a necessary part of the analysis, a solution was presented 
which developed the mode two and mode three stress intensity factors 
for an embedded elliptical crack subjected to nonuniform shear 
stresses. Additional equations were derived for the residual stress 
components on the crack surface which resulted from the freeing of the 
surfaces of the plate of stresses during the alternating technique. 

Results were computed and presented graphically for the problem in 
which a uniform shear stress was applied to the surface of the crack 
and parallel to the semimajor axis of the ellipse, while no initial 
loading was applied in the direction parallel to the semiminor axis. 
Crack shapes wore varied from a/2c ratios of 0.05 to 0.40 whereas 
crack depth ratios a/t were varied from 0.2 to 0.9. When compared to 
the maximum K, and K. values for embedded ellipses, the results 
revealed that there was only a slight change in the maximum K 2 and Kj 
values for similar semielliptical cracks having low a/t ratios, 
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both and became substantially larger as the crack depth a/'t 
was increased, although experienced a significantly greater increase 
than did The maximum value of the mode three stress intensity 

for a crack with a/2c of 0.2 increas'd by 15.1 and 31.2 percent as a/t 
was increased to 0.8 and 0.0, respectively. The maximum mode two 
stress intensity factor for the corresponding geometries increased 
by 6.0 and 7 .\ percent, respectively. Similar increases occur for 
other problem geometries, so these factors can play an important role 
in situations where mixed mode loading can occur. 

The normal loading problem was solved in a manner which was 
essentially identical to the one used for the shear problem. The 
alternating method was again used to obtain the stress intensity 
factors, although it was the mode one factor Kj which was of interest 
in this case. 

Computer calculations were made and results were plotted for the 
case in which 3 uniform normal stress was applied to the surface of 
the crack in a direction which was perpendicular to the plane of the 
crack. Emphasis was placed on slender cracks, so runs were made for 
cases in which the crack shapes a/2c varied from 0.05 to 0.3 as crack 
depth ratios a/t were varied from 0.2 to 0.95. 

The maximum value for the mode one stress intensity factor 
increases for any particular crack depth a/t as the crack becomes more 
slender with a decreasing a/2c ratio. Likewise, K increased very 
significantly a*, tlu* a/t ratio was increased for any particular a /2c 
value. Very slender cracks of low n/2c ratios and high a/t ratios 
produce high mode one stress intensity factors, and so are of special 
importance to the designer. Additionally, the maximum stress intensity 
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does not always occur at the deepest part oi the crack, but can occur j 

i 

up to 40 or 50 degrees from that point. 

Results from this study were compared with those obtained by other 

researchers who neglected the effects of the interaction between the | 

i 

front and back surfaces. The results compare well for cracks having ! 

i 

a/2c of 0.2 and greater, but some difference in the results of the two j 

> studies is apparent for lower a/2c ratios, indicating the importance 

of this effect. 

/ 

The results of this study were used to calculate the fracture j 

toughness values for a large number of epoxy plates on which mode one j 

tests were made. When normalized, these fracture toughness values are ! 

\ 

nearly constant, thus indicating the validity o r the calculated results. i 

i 

2.0 Recommendations 

i m ■ ■ — - I 

The complex loading of today’s structures makes either of the 
two loadings discussed incomoletc or inadequate when considered alone. 

i 

When more than one mode is acting along the crack border, failure of 
the part becomes a function of all nonzero stress intensity factors. 

Very little work has been done on developing a theoretically or experi- 
mentally obtained mixed mode fracture criterion. This criterion is an 
urgent requirement if information such that obtained xn this study 
is to be employed usefully in the field of applied fracture mechanics. 

# As presently developed, the alternating method is limited to 

* geometric configurations in which the front and back surfaces are 

parallel. This restriction on the relative placement of the surfaces 
does not allow the solution to be obtained for a large class of 
problems. A new techniiu e utilizing a better method of obtaining the 
solution to stresses in an unerackcd body of general shape with general 
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boundary conditions, >uch as a three-dimensional finite element ^thod, 
would greatly incrct »e the problem solving capabilities of this method 
if it could be done without extensive increases in computer time 
requirements. 


* 
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Figure 9. Stress Intensity Facto 
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Figure 10. .‘tress Intensity Factors for a Surface Crack with Uniform 
Shear Applied Parallel to the Major Axis; a/2c * 0.30 
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Stress Intensity Factors for a Surface Crack with Uniform 
Shear Applied Parallel to the Major Axis; a/2c * 0.20 
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Figure 12. Stress Intensity Factors for a Surface Crack with Uniform 
Shear Applied Parallel to the Major Axis; a/2c = 0.10 
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figure 1(>. Stress Intensity factors for Successive Cycles of 
Iteration; n/2c = 0.20 and a/t ~ 0.90 
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Artungcmcnt of Surface Rectangles for the Normal 
Loading Problem 
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Figure 18. Stress Intensity Factor for a Surface Crack Subjected 
to Normal Loading; a/2c * 0.30 
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Figure -0. Stress Intensity Factor for a Surface Crack Subjected 
to Normal loading; a/2c * 0.10 
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Figure 2 &. Fracture Toughness Data for Epoxy 
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APPENDIX I 


The non zero terms of Liu* matrix [K^] are ^Wen by the expression 
below and the parameters A, B, C, p and q are given in the table. 


K. , « A[B(l-v)M + CvM ] 2G 

ij p q 


















































p 

a 

P 

1. 

2/a 3 k 2 

2. 

2/ab 2 k 2 

3. 

2/aV 2 

A. 

2/a Vk 4 

5. 

2/a 5 k 2 k' 2 

6. 

2/a s k z k' 4 

7. 

2/a 7 k 4 k' 4 

8. 

2/3a 5 b 2 k 6 

9. 

2/3a 3 b 4 k 6 

10. 

2/3a 7 k 4 k' 2 

11. 

2/3a 7 k 4 k' 1 

12. 

2/3a 7 k 2 k* 4 

13. 

2/3a 7 k 2 k' 6 

14. 

2/3a s k 4 

15. 

2/3ab 4 k 4 

16. 

2/15a 7 k 6 

17. 

2/15ab 6 k 6 

18. 

2/15a 7 b 2 k 8 

19. 

2/3a 5 b 4 k e 

20. 

2/15a 3 b 6 k 8 

21. 

2/15a 9 k 6 k' 2 

22. 

2/15aVk’ 8 

23. 

2/3aVk’ 6 
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APPENDIX I 
(continued) 


M p - ct p [6 p E(k) + Y p K(k) ] 


8 p 

Y P 

-1 

+1 

+1 

-k' 2 

-1 

0 

1+k 1 2 

-2k' 2 

l-2k"’ 

-k’ 2 

-1-k 2 

k' 2 

-2(k’ 2 +k 4 ) 

k' 2 (2-k 2 ) 

3+5k ,2 +2k 2 k ,z 

-k' 2 (8+k 2 ) 

3k 2 k' 2 +10k 2 -8 

8k ,2 -9k 2 k' 2 

2+3k 2 -8k 4 

-2(2k 2 +l)k' 2 

3k 2 +(2-10k 2 )/k' 2 

2 ( 3k. 2 —1 

8k 2 k' 2 +5k 2 -3 

(3k' 2 -k 2 )k ,z 

2k 2 k ,2 +5k 2 +3 

~(k 2 +3)k' 2 

-2(l+k 2 ) 

?+k 2 

2 ( 2k 2 — 1 ) 

k' 2 (2-3k 2 ) 

-(8k 4 +7k 2 +8) 

4k 4 +3k 2 +8 

23< 4 -23k 2 +8 

-k ,2 (15k 4 -19k 2 +8) 

15+33k’ 2 +.L7k 2 k’ 2 +8kV 2 

-Ak' 2 (k 4 +2k z +12) 

18k 2 -16+8k 2 k ' ? -2k 2 k' 4 

k' 2 (16-17k 2 +k 2 k’ 2 ) 

73k 4 -113k 2 +48-15k 2 k' 4 

-k' 1 (60k 4 -104k 2 +48) 

8+9k 2 +16k 4 -48k 6 

-k' 2 (2Ak 4 +13k 2 +8) 


40k 2 k’ 2 -48k 2 +88k V 2 -15k 4 k' 4 k' 2 (16k 4 -45k V 2 +8k 2 ) 

llk 2 k' 2 -8k 2 -k 2 k' 4 


2k 2 k' 4 +Ak 2 k' 2 -24k 2 -16+16/k' 2 
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APPENDIX I 
(continued) 


p 

a 

P 

6 p 

Y P 

24. 

2/3aVk , “ 

16k 4 k ,2 +8k 4 -4k 2 ~2 

k' 2 (2+5k 2 -8k 4 ) 

25. 

2/105ab 8 k e 

8(2k 2 -l)(6-llk 2 k' 2 ) 

k’ 2 (105k 4 k' 2 +88k"-160k 2 +48) 

26. 

2/3a 9 kV 6 

8-2k ,2 -3k 2 k ,z -8k 2 k’ 4 

2k* 2 (2k z k' z -k f z -2) 

27. 

2/3a 9 k 6 k' 4 

8k 4 k ,2 +3k 2 k ,2 -2k 2 -8 

k' 2 (4k z k* 2j -5k’ 2 +3) 

28. 

2/3a 9 k 6 k' 4 

6k 4 -3k 2 +16-8/k ,z 

3k 2 k ,z +12k 2 -8 

29. 

2/105a 9 k 8 

8(l+k 2 )(-6k 4 +k 2 -6) 

24k 6 +I7k 4 +16k 2 +48 
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APPENDIX II 


If 


The following are the expressions for K and K corresponding to 

2 3 

each term iu the stress function series. 


a and b 
0 0 0 0 

K = 1 7 r- (a b cos(J) + b a sin<J>) 

2 (ab) 3/2 0 0 0 0 


K s QQH (1-v) (b b cos<j> - a a sin$) 

3 (ab) 3/z <T 00 00 

a and b 
to 10 

% 

K m ^2G)t ( a t C os 2 4> + b sin4> cos<J>) 

2 (ab) 3 / ? ** 10 a 10 

K a 3? — w — (i-v) (a alnij) costj) - b — cob 2 <Ji) 

3 (ab) 3 '' 2 or 10 10 a 


a and b 
01 01 

% 

K « — ( a sin<|> cos<J> + b 7 ~ sin 2 <J>) 

2 (ab ) 3 / 2 01 01 b 


K - ? - 9 ~- (l-v) (a z sin 2 <f> - b sinifi cos(|)) 

3 (ab) 3 / 2 or 01 b 01 


a and b 

11 1 1 

2 b 

K ■ - r- ( — - sini|> cos ? i(i + r — sin'’<|> eo8<j>) 

2 (ab) 3/2 A 8 b 


K - _ (i-v) ( — ~ 3in4> co 8 2 4) - r— sln 2 4- cob4>) 

3 (ab) 3/ 2 0 T a b 


■ 





a and b 
? o 20 

„ , H2GH 1 * r ( , S 2 cos'* + b - sin* cos 2 *) 

‘ 2 (ab) 2 ' 2 * 11 2 " a 

k - — c a-v) (b 20 **«>•’♦ - ‘J sln * cofl! * ) 

3 (ab) 3/2 * 


a and b 

02 02 

I92Gtt__ (a I C os<fr + b “2 sin 3 <}>) 

K “ IT, a ^02 b ^ 02 b 

2 (ab) 3/2 ♦ 

K . ±M]£ (i-v) (b £ sin 2 <(> cos<|> - a Q2 p sln3<,>) 

3 (ab) 3 / 2 r 02 


a and b 

2 1 21 

. ;J53bClI* I 2 .in* cos’* + b JJ 8lnJ * C ° s!,H 

2 (ab)’/’* 14 ” “ 

. 15 . 360^ (1 . v) (a ij aln 2 * co. 2 * - b^ ±2 »i»* 

> (ab)* /2 ♦' 


a and b 

12 12 

-lSSbCb 12 _ u L sln 2 * co. 2 t + b J! »in’4> 

2 (ab)*' 2 ** 12 * b 

K . 133bGoj_ (1 . v , ( i 2 .in’* co.* - b^ ^ -l- 2 * co. 2 ») 

3 (ab ) 3/2 ^ W 


3 and b 

30 30 

. -1536G1TL. b co.'* + b Ji >tn* »■*♦> 

2 (ab)*' 2 * 11 “ 

. 153 W 1 * ( 1 -y) (a 4 .lb* CO.’* - b J> co.'*) 

* (ab )*' 2 ** 


\)7 


a and b 

0 3 0 3 


K * — r (a t “2 sin 3 <j) cos<}> + b r3 sin 4 (J>) 

2 (ab) 3 ' 2 P 03 b 03 


K - 


1536GTT' 


3 (ab) 3/2 ^ 


(1-v) (a ^3 sin 4 <(i - b h sln 3 <p cos<J>) 
0 3b 03» 


4> * a 2 sln 2 <p + b 2 cos z (p 



